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FOREWORD 
(Nontechnical  summary) 

This  report  describes  a  numerical  integration  program  much  like  the  amended 
form  of  NUMIT-ONE  but  written  for  the  SDS  920,  a  faster  and  larger  computer. 
While  the  flow  chart  is  identical  to  the  amended  version  of  NUMIT-ONE,  the  coding 
has  taken  advantage  of  the  more  flexible  SDS  FORTRAN  II  and  the  input  and  output 
devices  have  been  changed  to  paper  tape  since  this  is  standard  equipment  for  SDS  900 
series  computers  and  a  card  read -punch  is  not  available  for  our  SDS  920. 
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ABSTRACT 

An  amendment  to  NUMIT-ONE  is  presented.  The  program  utilizes  Simpson's 
Rule  to  perform  integration  over  any  closed  intervals  on  the  SDS  920. 
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I.  INTRODUCTION 


NUMIT-ONE,  originally  coded  for  the  IBM  1620,  has  now  been  recoded  for  the 
SDS  920.  The  code  follows  the  flow  chart  given  in  the  NUMIT-ONE  amendment.  * 
The  code  itself  approximately  duplicates  the  IBM  1620  code;  the  only  changes  being 
those  to  take  advantage  of  SDS  FORTRAN  II  and  input/output  instructions.  Familiarity 
with  both  of  the  previous  publications  is  assumed. 

Equipment  needed: 

FORTRAN  II  Processor 

SDS  920  with  minimum  4  K  memory 

While  it  has  not  been  attempted,  this  program  could  probably  be  compiled  and 
run  on  a  smaller  SDS  computer  of  the  900  series. 

II.  ACCURACY 

Floatingpoint  quantities  in  the  920  are  held  to  approximately  12  digits  unlike 
the  IBM  1620  where  variable  precision  parameters  are  possible.  Generally,  one  can 
obtain  six  significant  figures,  the  maximum  output  precision  now  called  for. 

Proper  analysis  prior  to  attempting  to  integrate,  e.g. ,  plotting  data  to  visually 
select  subintervals  of  integration,  will  have  the  double  benefit  of  improving  accuracy 
and  reducing  machine  time.  The  test  case  in  Appendix  B  illustrates  this.  On  the 
first  rim  the  total  interval  was  broken  into  three  parts ,  with  each  subinterval  carefully 
chosen  so  that  unnecessary  extremes  would  not  be  reached  in  any  interval  and  the 
results  for  the  parts  were  summed  after  the  three  passes.  The  result  is  an  answer 
accurate  to  approximately  six  significant  figures  in  a  few  minutes  of  machine  time. 
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One  can  integrate  over  the  entire  interval  in  cases  when  time  for  analysis  is 
short;  that  is,  when  plotting  and  other  manual  analysis  is  not  possible.  However,  the 
second  rim  of  the  test  case  in  Appendix  B  using  the  entire  interval  method  took  over 
twice  as  much  machine  time  and  the  result  was  slightly  less  accurate  than  the  previ¬ 
ous  run  after  14  iterations.  Analysis  is  important  when  time  permits. 

A  word  of  caution  about  the  limits  of  integration:  a  situation  in  which  the  pro¬ 
gram  will  break  down  occurs  when  integrating  over  an  interval  where  the  integrand 
takes  on  both  positive  and  negative  values  and  produces  an  integral  close  to  zero.  In 
such  cases ,  significant  figures  will  be  lost  during  integrating  and  DELTA  J  will  exhibit 
a  wild  behavior  instead  of  converging  rapidly  as  it  should.  In  a  case  like  this,  it  is 
best  to  split  the  integrand  into  at  least  two  parts,  each  containing  a  nonzero  part  of 
the  total  integral. 

III.  INPUT  LIMITATIONS 

Input  limitations  are  the  same  as  for  the  amended  NUMIT-ONE  for  the  IBM  1620 
with  one  exception.  Jmax  f°r  this  version  must  not  be  more  than  24  due  to  limits  on 
the  size  of  fixed  point  quantities  in  SDS  FORTRAN  II.  Generally  this  allows  for  more 
iterations  than  one  would  be  willing  to  use. 

The  integrand,  f(x),  is  provided  by  a  subprogram  of  the  FUNCTION  type  called 
FUN  (X,N).  X,  computed  in  the  main  program,  is  the  value  of  X  when  the  function  is 
called  on;  N*  is  the  "pass  number".  FUN  must  appear  on  the  left  of  an  arithmetic 
statement  equivalent  to  the  integrand  in  the  subprogram  before  returning  to  the  main 
program. 

*  The  integral  is  evaluated  over  several  separate  X  regions  as  a  sum  of  several  inde¬ 
pendent  integrals  over  disjoint  regions;  the  integer  variable  describing  the  separate 
integrals  is  N. 
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IV.  INPUT  AND  OUTPUT  FORMATS 


Input  Format 

Input  is  by  paper  tape.  The  following  quantities  are  specified  in  the  input  (float¬ 
ing  point  quantities  are  in  lower  case,  fixed  point  quantities  in  upper  case), 
a  =  lower  limit  of  interval 
b  =  upper  limit  of  interval 
Jmax  =  maximum  number  of  iterations  desired 
A  =  convergence  criterion  (floating  point) 

Line  1 

Alpha  Identification.  Reproduced  in  the  output  and  does  not  affect  process¬ 
ing.  Up  to  80  characters  permitted. 

Line  2 

Integrand  Specification.  Contains  alphanumeric  information  to  identify 
integrand  in  output.  This  does  not  affect  processing.  Up  to  48  characters  permitted. 


Line  3 

Columns  1-11 

a 

Format 

E  11.6 

Columns  12-22 

b 

E 

11.6 

Columns  23-25 

Jmax 

I 

3 

Columns  26-33 

A 

F 

8.6 

Columns  34-35 

Switch 

F 

2.0 

Following  line  3  may  be  any  number  of  additional  lines  with  new  values  of  a,  b, 
Jmax»  and  A  punched  according  to  the  same  formats.  If  more  lines  follow,  switch 
must  be  punched.  If  no  lines  follow,  switch  should  not  be  punched.  In  the  last  line, 
switch  must  not  be  punched. 
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Errors  in  Input 


All  the  error  checks  listed  in  AFRRI  SR65-1  are  made,  with  one  change.  Jmax 
is  permitted  to  be  as  high  as  24  instead  of  11  before  an  error  message  will  be  punched 
and  the  computer  halts . 

Output  Format 

Output  is  identical  to  the  amended  form  of  NUMIT-ONE  for  the  IBM  1620,  the 
only  change  being  that  output  is  on  paper  tape  rather  than  cards . 

Two  sample  listings  are  in  Appendix  B. 

V.  OPERATION 

Operation  is  routine.  The  standard  FORTRAN  II  compiler  is  used  and  no 
switches  are  utilized. 
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NUM1T  ONE 
D  IMENS  ION 

PC0]»1. 

P [1 ]  =  4  , 

P  E2]*l , 


-  ROCKWELL  -  NCV,  1965  , 

FC0/2],PC0/2],Q[0/2],U[6] 


Q  [0 ]-l * 

Q[1  ]  =  2  , 

0.  C23®1  • 

N=0 

TOTAL  ■  0, 

ACCEPT  TAPE  100 
PUNCH  TAPE  100 
PUNCH  TAPE  101 

ACCEPT  TAPE  107, CU [| ] , t *1 ,6] 

PUNCH  TAPE  108,  EU  E  l  ]  , 1=1/6] 

1  ACCEPT  TAPE  1 02 ,A,B,JHAX ,DELTA1  ,SW1 
I  Ft  I FtA,B] 33 ,2,3 

2  PUNCH  TAPE  103 
STOP 

3  |  F  EDE  LTAl  ]  5 , 4 ,5 

4  PUNCH  TAPE  104 


STOP 

5  | F  EJMAX -24]7,7,6 

6  PUNCH  TAPE  105 


STOP 

7  I F  EJMAX -2] 6,8,8 

8  SUMl=0, 

SUK-0. 

N*N+1 


PUNCH  TAPE  1 0  6,B , JMAX , A,DE LTAl 
9  DELX*  EB-A]/ E2 ] 

I  F  E  J-l  ]10  /I  0/1 4 

10  DO  11  Ic0,2 
X=A+  l*DELX 
FEU  =  FUN  EX,N] 

SUMl=SUMl+QE 1 ]*F  E I ] 

11  SUM  * S UM  +  P[|]*F[|] 
ANl=DELX*SUM/3% 

J  ■  J+l 
GO  TO  9 

14  SUM=0  t 

DO  15  |*1, K2 
X=A+  E2* I -1 ]*DELX 

15  SUM  «=  SUM+FUN  EX,N] 

AN2  *  DELX*E4.*SUM+SUMl]/3, 

DELTA  =  ABSFE EANl-AN2]/AN2*100t] 
I F  EDELTA-DELTAl ]1 9 ,19,16 

16  I FEJ-JMAX117 ,19,19 

17  ANl=AN2 

SUMl=  $UM1+2,*SUM 

PUNCH  TAPE  109, J,AN2, DELTA 

J*J+1 

K2C2*K2 
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«  56  GO  TO  9 

=  57  19  PUNCH  TAPE  110,J,N,AN2,DELTA 

»  5  8  T  CTALeT  OTAL+AN2 

-  59  I F[SW1]1 , 20,1 

=  60  20  PUNCH  TAPE  lll,N,TOTAL 

=  61  STOP 

=  62  100  FORMAT [80H 

*  63  1  ] 

=  64  107  FOPMATC6A8] 

=  65  108  FORMAT [15X, $F [X]  =  S,6A8] 

=  66  101  FORMAT  [/15X,  $NUf1  |T  ONE  [ONE-DIM.  INTEGRATION,  SIMPSONS  RULE]$] 

=  67  102  FORMAT [2E11, 6, |3 ,F8 .6 ,F2 .0 ] 

-  68  103  FORMAT [10X,$CHECK  INPUT  -  A«D$1 

=  69  104  F ORMAT [1 0X , $CHE CK  INPUT  -  DELTA  *  0ft ] 

■  70  105  FORMAT [1 0X , $CHE  CK  INPUT  -  JMAX  OUTSIDE  RANGE  $  ] 

=  71  106  FORMAT  [//  3  5X  „  ft  INPUT  DATA$/21X,$B  =  S1PE1 2  .5 , 9X,  $J  [MAX  ]  ■=  $13/ 

=  72  1  21X,  $A  =  $1PE1 2 ,5 , 9X  ,  SDE  LTAl  =  $1PE10 .4//1 6X,  $  |TE  RAT  I  ON  $,11X, 

»  73  2$|NTEGRAL$,13X,$DELTA$] 

»  74  109  F ORMAT [18X, |3,12X,1PE12,5/8X/0PF10.4] 

=  75  110  FORMAT [18X, 13 , $ [PASS$,  |2 , $]$ ,4X ,1PE1 2,5,8X/0PF10, 4/// ] 

=  76  111  FORMAT  [3 2X$NO,  PASSES  =  $  13/ 29X$TCTAL  INTEGRAL*  $1PE12.5/ 

=  77  139X$END$] 

=  78  END 


PROGRAM  ALLOCATION 


00007  F 

00015  P 

00023  Q 

00031  U 

00045  N 

00046  1 

00047  JMAX 

00050  J 

00051  K2 

00052  TOTAL 

00054  A 

00056  B 

00060  DE LTAl 

00062  SW1 

00064  SUM1 

00066  SUM 

00070  DELX 

00072  X 

00074  .AN1 

00076  AN2 

00100  DELTA 

SUBPROGRAMS  REQUIRED 


IF  FUN  ABSF 

THE  END 
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APPENDIX  B 


Test  Case 

Running  time  is  still  highly  variable  and  depends  strongly  on  the  integrand. 
However,  time  is  much  shorter  than  on  the  IBM  1620.  No  case  is  known  which  is  not 
at  least  15  times  as  fast  on  the  SDS  920  as  on  the  other  computer,  the  ratio  being 
highest  when  the  number  of  iterations  is  large.  The  formula  given  in  the  previous 
reports , 

t  =  K  •  2Jact 


where  K  is  a  constant  dependent  on  the  integrand,  still  appears  to  be  valid  when  the 
number  of  iterations  is  large  enough  so  the  time  spent  computing  is  large  compared 
to  the  time  spent  punching  tape. 

The  function  illustrated  here, 

m 


=  y7  —  l-X2 

(2-x)13/2 


has  been  the  most  difficult  tested  thus  far.  It  converges  slower  and  the  value  of 
K  (7.47xl0-4sec)  is  larger  than  any  heretofore  reported.  Figure  1  illustrates  the 
function  over  the  intervals  of  -1  to  1.  Running  time  for  RUN  1  was  5.5  minutes; 
RUN  2  took  12.25  minutes. 


9 


Figure  1 
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C  SOURCE  STATEMENTS  OF  SUBROUTINE 
C 

FUNCTION  FUN(X.N) 

FUN=X** 7*SQRT  F (1 . -X*  X )/ ( 2 . -X ) **  6 . 5 

RETURN 

END 


INPUT  DATA 
RUN  1 


,  ,  ,  TEST  FUNCTION 

X**7*SQRT(l-X*x)/(2-X)**13/2 
- 100000+01, +500000+00, 12,. 010000,1, 
+500000+00, +800000+00, 12, .001000, 1 , 
+800000+00, +100000+01 ,12,. 000100, 


RUN  2 


X**7*SQRT( 

-100000+01,+ 


1-X*X)/(2-X)**l3/2 
100000+01,14,. 000500 


TEST  FUNCTION 
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OUTPUT  DATA 


RUN  1 

TEST  FUNCTION 

NUMIT  ONE  [ONE-DIM.  INTEGRATION,  SIMPSONS  RULE] 
F[X]  =  X**7*SQRT(l-X*X)/(2-X)**13/2 


INPUT  DATA 

B  =  5.00000E-01  J[MAX]  =  12 

A  =  -1.00000E  00  DELTA1  =  9.9999E-03 


ITERATION 

INTEGRAL 

DELTA 

2 

3.31206E-05 

265.1566 

3 

-1.74285E-05 

290.0372 

4 

-3.03975E-05 

42.6647 

5 

-3.29935E-05 

7.8681 

6 

-3.36745E-05 

2.0223 

7 

-3.38945E-05 

0.6491 

8 

-3.39703E-05 

•0.2231 

9 

-3.39968E-05 

0.0781 

10 

-3.40062E-05 

0.0275 

11 [PASS  1] 

-3.40095E-05 

0.0097 

INPUT  DATA 

B  =  7.99999E-01  J[MAX] 

A  =  5.00000E-01  DELTA1 


12 

9.9999E-04 


ITERATION  INTEGRAL  DELTA 

2  2.90366E-03  3.5586 

3  2.89730E-03  0.2193 

4  2.89694E-03  0.0126 

5 [PASS  2]  2.89692E-03  0.0008 


INPUT  DATA 

B  =  1.00000E  00  J[MAX]  =  12 

A  =  7.99999E-01  DELTA1  =  9.9999E-05 


ITERATION 

2 

3 

4 

5 

6 

7 

8 
9 

10 
1 1 

12[PASS  3] 


INTEGRAL 
1 .95064E-02 
2.0507  IE-02 
2.08279E-02 
2.09361E-02 
2.09735E-02 
2.09866E-02 
2.099 12E-02 
2.09928E-02 
2.09934E-02 
2.09936E-02 
2.09936E-02 


DELTA 
16.7288 
4.8794 
1.5404 
0.5167 
0.1784 
0.0624 
0.0219 
0.0077 
0.0027 
0  .00  1 0 
0.0003 


NO.  PASSES  = 
TOTAL  INTEGRAL= 
END 


3 

2.38566E-02 
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RUN  2 


TEST  FUNCTION 

NUMIT  ONE  [ONE-DIM,  INTEGRATION,  SIMPSONS  RULE] 
F[X]  =  X**7*SQRT(l-X*X)/(2-X)**13/2 


B  =  1 

INPUT  DATA 
.00000E  00 

J[MAX]  =  14 

A  =  -1 

.00000E  00 

DELTA  1  =  4.9999E-04 

ITERATION 

INTEGRAL 

DELTA 

2 

3.11637E-04 

100.0000 

3 

6.93780E-03 

95.5081 

4 

1.70652E-02 

59.3454 

5 

2.17058E-02 

21.3795 

6 

2.31639E-02 

6.2948 

7 

2.36230E-02 

1.9433 

8 

2.37758E-02 

0.6426 

9 

2.38283E-02 

0.2204 

10 

2.38466E-02 

0.0769 

1 1 

2.3853  IE-02 

0.0270 

12 

2.38554E-02 

0.0095 

13 

2.38562E-02 

0.0034 

14[PASS 

1]  2.38564E-02 

0.0012 

NO.  PASSES  =  1 

TOTAL  INTEGRAL=  2.38564E-02 
END 
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